load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"

begin
  YYYY  = ispan(1951,2100,1)
  tStrt = ind(YYYY.eq.1995)
  tLast = ind(YYYY.eq.2014)

  model = (/"ACCESS-ESM1-5","MIROC6","ACCESS-CM2","BCC-CSM2-MR","CNRM-CM6-1","GFDL-ESM4",\
             "HadGEM3-GC31-LL","IPSL-CM6A-LR","MRI-ESM2-0","NorESM2-LM"/) 

  pr    = new((/10,150/),float)

  do m=0,9
    fread = addfile("data/pr_timeseries_1951-2020_JJAS_historical_"+model(m)+".nc","r")
    preh  = fread->prN; averaged precipitation over northern HMA

    fread = addfile("data/pr_timeseries_2021-2100_JJAS_ssp245_"+model(m)+".nc","r")
    pres  = fread->prN

    n  = dimsizes(preh)
    pr = new((/n(0),150/),"float")

    pr(:,:69) = preh 
    pr(:,70:) = pres

    pr   = pr-avg(pr(:,tStrt:tLast))
    pr   = runave_n_Wrap(pr,9,0,1)    

    t    = ispan(1,142,1)
    do i=0,n(0)-1 
      cc = lspoly(t,pr(i,4:145),1,5); Missing values are excluded
      pr(i,4:145) = cc(0)+cc(1)*t+cc(2)*t^2+cc(3)*t^3+cc(4)*t^4
    end do 

    pr(m,:) = dim_avg_n_Wrap(pr,0); The signal for each model

    delete([/pr,preh,pres/]) 
  end do 

end 